/*
 * Gurobi.cpp
 *
 *  Created on: May 14, 2013
 *      Author: sjelic
 */

#include "Gurobi.h"


namespace coveringPacking {


	int __stdcall mycallback(GRBmodel *model, void *cbdata, int where, void *usrdata)
	{
		struct testData *mydata = (struct testData *) usrdata;
		if (where==GRB_CB_BARRIER)
		{
			GRBcbget(cbdata, where, GRB_CB_BARRIER_PRIMOBJ,&mydata->primobj);
			GRBcbget(cbdata, where, GRB_CB_BARRIER_DUALOBJ,&mydata->dualobj);

			if(mydata->primobj!=0)
			{
				if((mydata->dualobj/mydata->primobj)<1+mydata->eps&&(mydata->dualobj/mydata->primobj)>=1)
				{
					GRBterminate(model);
				}
			}
			else
			{
				if(mydata->dualobj<mydata->eps&&mydata->dualobj>=0)
				{
					GRBterminate(model);
				}
			}


		}
		return 0;
	}
	template <typename T>
	T** AllocateDynamicArray( int nRows, int nCols)
	{
		  T **dynamicArray;

		  dynamicArray = new T*[nRows];
		  for( int i = 0 ; i < nRows ; i++ )
		  {dynamicArray[i] = new T [nCols]();}

		  return dynamicArray;
	}
	template <typename T>
	T** ZeroOutDynamicArray( int nRows, int nCols, T** A)
	{
		for(int i=0; i<nRows; i++)
		{
			for(int j=0; j<nCols; j++)
			{
				A[i][j]=0;
			}
		}
		return A;
	}

	template <typename T>
	void FreeDynamicArray(T** dArray)
	{
		  delete [] *dArray;
		  delete [] dArray;
	}

	GurobiLP::GurobiLP(double D, double EPS, string optimizingLogFile, bool NORMALIZED)
	{
		GRBenv* env = NULL;
		GRBmodel* myModel = NULL;
		myTestData.m=0;
		myTestData.n=0;
		myTestData.d=D;
		myTestData.eps=EPS;
		myTestData.sep='&';
		logFilePath=optimizingLogFile;
		normalized=NORMALIZED;
	}
	void GurobiLP::convertMatrixAFromYoungToGurobiModel(string path)
	{
		//-----------------------------------------------------------------------------------
		FILE *filePt;
		filePt = std::fopen(path.c_str(),"r");
		unsigned int M, N, numberOfNonzeroElements;

		if(!filePt){ printf("Input file missing ...\n"); exit(0);}
		if(std::fscanf(filePt,"%d %d %d",&M,&N,&numberOfNonzeroElements)){};
		myTestData.m=M;
		myTestData.n=N;
		myTestData.numOfNonZeroInConstrMtrx=numberOfNonzeroElements;


		matrix_type *A;
		A = (matrix_type *) calloc(M*N,sizeof(matrix_type));
		index row=0,col=0;
		integer readings=0;
		while(readings++<numberOfNonzeroElements)
	   {
			matrix_type value;
			if(std::fscanf(filePt,"%d %d %f",&row,&col,&value)){};
			A[col+row*N] = (matrix_type) abs(value); // row-wise storing
	   }


		if(!normalized)
		{
			real *b = (real *) calloc(M,sizeof(real));
			real *w = (real *) calloc(N,sizeof(real));

			index row = 0,col = 0;
			integer readings = 0;

			// read array c
			readings=0;
			while(readings<N)
			{
				float value;
				if(std::fscanf(filePt,"%f",&value)){};
				w[readings] = (real) value;
				if(readings==0) std::printf("%.3f",value);
						readings++;
			}
			// read array b
			readings=0;
			while(readings<M)
			{
				float value;
				if(std::fscanf(filePt,"%f",&value)){};
				b[readings] = (real) value;
				if(readings==0) std::printf("%.3f",value);
				readings++;
			}


			for(int i=0; i<M; i++)
			{
				for(int j=0; j<N; j++)
				{
					A[i*N+j]/=w[j]*b[i];
				}
			}
			free(b);
			free(w);


		}
		std::fclose(filePt);




		//-----------------------------------------------------------------------------------

		int *vbeg =  new int[N](); 		//memory allocation for sparse representation of constraint matrix in Gurobi C
		int *vlen = new int[N]();
		int *vind = new int[numberOfNonzeroElements];
		double *vvalue = new double[numberOfNonzeroElements];
		char *sense = new char[M];
		double *rhs = new double[M];
		double *obj = new double[N];
		char *vtype = new char[N];
		for(int j=0; j<N; j++) obj[j]=1;//initialization
		for(int j=0; j<N; j++) vtype[j]=GRB_CONTINUOUS;
		for(int j=0; j<N; j++) vbeg[j]--;
		for(int i=0; i<M; i++) sense[i]=GRB_LESS_EQUAL;
		for(int i=0; i<M; i++) rhs[i]=1;

		int count=0;
		for(int j=0; j<N; j++) //calculation of sparse representation of constraint matrix
		{
			for(int i=0; i<M; i++)
			{
				if(A[i*N+j]!=0)
				{
					if(vbeg[j]==-1){vbeg[j]=count;}
					vlen[j]++;
					vind[count]=i;
					vvalue[count]=A[i*N+j];
					count++;
				}
			}
		}

		//-----------------------------------------------------------------------------------

		if(GRBloadenv(&env,logFilePath.c_str())) std::cout<<"Loading environment error!"<<std::endl; //loading environment, model and updating model
		if(GRBloadmodel(env, &myModel, "Covering Packing Example", N, M, -1, 0.0, obj, sense, rhs, vbeg, vlen, vind, vvalue, NULL, NULL, vtype, NULL, NULL)) std::cout<<"Problem with loading Gurobi model"<<std::endl;
		if(GRBupdatemodel(myModel)) std::cout<<"Problem with updating model"<<std::endl;

		delete[]obj; //deallocation of memory
		delete[]sense;
		delete[]rhs;
		delete[]vbeg;
		delete[]vlen;
		delete[]vind;
		delete[]vvalue;
		delete[]vtype;
		free(A);
	}

	void GurobiLP::tuningAndSolvingModel()
	{
		//if(GRBsetdblparam(GRBgetenv(myModel),"BarConvTol",myTestData.eps)) std::cout<<"Problem with tuning BarConvTol parameter."<<std::endl;
		//if(GRBupdatemodel(myModel)) std::cout<<"Problem with updating model"<<std::endl;
		if(GRBsetdblparam(GRBgetenv(myModel),"FeasibilityTol",1e-9)) std::cout<<"Problem with tuning FeasibilityTol parameter."<<std::endl;
		if(GRBupdatemodel(myModel)) std::cout<<"Problem with updating model"<<std::endl;
		if(GRBsetintparam(GRBgetenv(myModel),"Method",2)) std::cout<<"Problem with tuning Method parameter."<<std::endl;
		if(GRBupdatemodel(myModel)) std::cout<<"Problem with updating model"<<std::endl;
		if(GRBsetintparam(GRBgetenv(myModel),"Crossover",0)) std::cout<<"Problem with tuning Crossover parameter."<<std::endl;
		if(GRBupdatemodel(myModel)) std::cout<<"Problem with updating model"<<std::endl;
		//if(GRBwrite(myModel,"/home/sjelic/coveringPackingLinearSolver/5000x5000-0.5.lp")) std::cout<<"Problem with writing model to file"<<std::endl;
		if(GRBoptimize(myModel)) std::cout<<"Problem with optimizing model"<<std::endl;
		int status;
		if(GRBgetintattr(myModel,"Status",&status)) std::cout<<"Problem with getting Status parameter."<<std::endl;
		if(status==2 && myTestData.primobj!=0)
		{
			myTestData.primDualGap=myTestData.dualobj/myTestData.primobj;
		}
		else
		{
			myTestData.primDualGap=0.0;
		}
		return;
	}

	void GurobiLP::writingResultsToFile(string pathToFileNameOfResults, char s)
	{
		if(s) myTestData.sep=s;
		int error;
		//void *cbdata;
		std::ofstream resultFile;
		resultFile.open(pathToFileNameOfResults.c_str(),std::ios::app);

		if(resultFile.is_open())
		{
			GRBgetdblattr(myModel,"Runtime",&myTestData.time);
			resultFile.width(7);
			resultFile.setf(std::ios_base::right);
			resultFile<<myTestData.m<<" "<<myTestData.sep;

			resultFile.width(7);
			//resultFile.setf(std::ios_base::right);
			resultFile<<myTestData.n<<" "<<myTestData.sep;

			resultFile.width(6);
			//resultFile.setf(std::ios_base::right);
			resultFile<<myTestData.d<<" "<<myTestData.sep;

			resultFile.unsetf(std::ios_base::right);
			resultFile.setf(std::ios_base::left);
			resultFile<<" ";
			resultFile.width(6);
			resultFile<<myTestData.eps<<myTestData.sep;

			resultFile.precision(5);
			resultFile<<" ";
			resultFile.width(13);
			resultFile<<myTestData.time<<myTestData.sep;

			resultFile<<" ";
			resultFile.width(13);
			resultFile<<myTestData.primobj<<myTestData.sep;

			resultFile<<" ";
			resultFile.width(13);
			resultFile<<myTestData.dualobj<<myTestData.sep;

			resultFile<<" ";
			resultFile.width(13);
			resultFile<<myTestData.primDualGap;
			resultFile<<std::endl;
		}
		resultFile.close();
		return;
	}

	GurobiLP::~GurobiLP()
	{
		GRBfreemodel(myModel);
	}

} /* namespace coveringPacking */
